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1. INTRODUCTION 


It is well known that the gray level of image's pixels 
depends on the radiance coming, in a particular spectral 
range, from the ground areas corresponding to the pixels. 
If each area is homogeneous no problem exists in its in- 
terpretation; if, on the contrary, a lot of different 
objects are contained, they contribute in different man- 
ners to the global gray level of the pixel. 

The problem becomes more important as greater dimension 
pixels are involved, i.e. the problem must to be faced 
to understand the informative content of the data. 

Images from HCMM are typical example of this with their 
pixels 500 X 500 m2 dimensioned. It is then necessary to 
study the so called scale effect to understand how deta^ 
led informations aggregate to give rise to low or very low 
resolution data as they are detected for example from a 
high altitude platform. 

The first question to which this work was devoted was the 
simulation of HCMM data starting from better resolution 
ones . 

The following approach was however adopted: 

- determine a methodology which permits the correlation 
of different resolution data; 

- test the methodology for the simulation of little al- 
titude increase (usage of aircraft data) ; 

- simulate if possible HCMM images. 

The methodology is based on filtering techniques in the 
hyphothesis that all changes due to altitude increase can 
be represented by a point spread function. These aspects, 
together vvrith a '*practical'’ introduction to the filtering 
technique, are particularly pointed out all along the 
work. As will be seen, no attempt is made however to si- 
mulate HCMM data; attention is mainly focused on the sea 
le effect in its general aspects also for its importance 
in problems such as registration between images and so on. 


2. SCALE EFFECT; THEORETICAL MODEL 
2»1. Possible reduction methods 

At first sight the simulation of an altitude increase 
of acquisition platform consists of an image's reduction 
and this is equivalent to a suitable smoothing. This can 
be achieved by two methods : 

a) neighhorhood averaging techniques 

b) filtering techniques 

a) Neighborhood averaging techniques. These consist, as 
it's well known, in an arithmetical or weighted average 
of the gray values on a number of pixels equal to reduc 
tion factor (for example 4x4 pixels if the reduction fac 
tor is 1G;1). 

b) Filtering techniques. If this operation is realised 
in the spatial domain, it consists of the convolution of 
the image with the weights mask and then of the decima- 
tion of the points according to the appropriate scale 
factor. 


It's interesting to note that the former approach can be 
viewed as a particular case of the filtering technique 
where the filter dimensions are equal to the reduction 
factor. In this case the t\v'o steps of convolution and decima 
tion give an output equal to the averaging approach out- 
put. 

In both cases the weights mask has the role of a point-spread 
function (RSF) (for a survey on PSF and MTF see appendix 1) 
so its shape (and hence the MTF shape) should represent the 
physical characteristics of the process. 

The averaging approach is certainly very convenient from 
a computer time point of view, but it is considerably 
limited for accurate simulations. The weights mask must 
have in fact a Fourier transform as similar as possible 
to the process itself MTF. Since the more PSF samples 
there are, the better the transform shape can be con- 
trolled, it is convenient to use large dimensioned masks. 

In averaging technique, however, mask dimensions are fi- 
xed by the scale factor so that for little reductions 
very little extent mask must be used with the consequent 
impossibility of controlling the MTF shape. 


2,2. Discrete image formation process 

To understand how the function simulating a scaling up 
can be obtained, it is convenient to examine the mat he 
matical processes generating discrete images. In a gene 
ral manner, the imaging process at the two altitudes "" 
can be briefly represented as in Fig . 2.1. 



Fig. 2.1 
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is the discrete image at lower altitude 
is the discrete image at higher altitude q, 
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is the global PSF for image formation at 
altitude q^^^ 

is the global PSF for image formation at 
altitude q^ * Kq^^ (K > 1) 


This means that the same ground object is seen in dif- 
ferent manners (i.e. with different PSFs) and so dif- 
ferent images are obtained. 


Fig. 2.2. shows a more detailed representation of discre- 
te image formation process. sampl 
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i.e. each PSF is seen as composed of a continuous PSF 
followed by a sampler. 

The assumption is made that the signal-to-noise ratio 
value is sufficiently high so that the noise can be igno 
red. 

It's interesting to carry out a graphical derivation (in 
one-dimensional form) of the process. 

The assumption is made that SIPSF are involved so that 
the degrading process can in the same way be represented 
by its MTF. 

The ground object corresponding to the image being con- 
sidered is represented by the function o(x) (Fig. 2.3) 
while it is assumed that the continuous PSF h-i of Fig. 
2.2 is given by Fig. 2.3b. The integration performed on 
the photodetector surface is mathematically I'epresented ^ 
by the convolution of o(x) with hi(x) (i.e. the multipn 
cation of 0(f) and Hjj^(f) (Fig. 2,3c) (x). 

The assumption can be made that the convolution output 
has its spatial frequency content limited to a fA freque 
cy. The filtered image is then spatially sampledv The as 
sumption can be made that the sampling system works so^ 
that the image has neither overlaps nor gaps between pi- 
xels. Referring, for example, to the scan in the flight 
direction, this condition is satisfied by the frequency 
[Ref, ij , 

f - 

na 

Where V = Platform velocity 

h = Platform altitude 
n = Number of sensors 


(x) It must be noted that the object function o(x) is 

obtained by means of truncation of an extended ground 
object function with a window function which is unity 
in the range T^ and zero elsev\rhere. This is equiva- 
lent to convolve the object transform with the win- 
1 

dow transform (i.e. a function and tliis intro- 

duces a certain amount of ripple; this is however, 
disregarded in this treatment). 


» 
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a » Istantaneous field of viewCIFOV) 

With this assumption and supposing the sampler is an ideal 
Dirac delta function, it can be written 

f^«P . I 6(X-j (1) 

^ S M 00 * 


where is the pixel dimension depending on the 

altitude and system IFOV. It is assumed that for q, al- 
titude T (See Fig. 2.3d). 

As it is stated in Fig. A, I. 2 the sampling is equivalent 
to a multiplication in space domain, i.e. a cconvolution 
in frequency domain. The result is shown in Fig. 2.3e. It 
must be noted that if the folding frequency (fg/2) is less 
than fA (i.e. if data are undersampled), the aliasing ph£ 
nomenon appears. 

To provide a mapping between transform domain and image 
space and to permit transformed data to be handled by com 
puter, it is necessary to discretize also Fourier domain. 
The so-called Discrete Fourier Transform (DFT) is accom- 
plished by sampling the continuous Fourier transform by 
means of a fsaii^P with 1/T frequency step (Fig. 2.3f); this 
produces the iteration of°the space signal as it appears 
from Fig. 2.3g. For the two-dimensional case, this perio 
dicity is equivalent to the situation depicted in Fig. 2.4. 



From this it follows fRef. 2, pg. 20j that when fini- 
te extent discrete arrays are convolved, "circular con- 
volution" must be used. 

It must be noted that in the frequency domain a number 
of samples is generated equal to spatial domain one (for 
primary components) . This number is obviously .M » 'T /.j, . 


The same graphical derival; ion can be carried out for up 
per altitude (q^) image, i, The input (ground object 
o(x)) is the saAe than in ’’Fig. 2.3a. Generally the acqu^ 
sition system is a different one so a new PSF, h 2 (x), is 
now involved; referring to the case of a loss of “resolu- 
tion, the width of its non- zero region is however greater 
than the lower altitude one . 


In the particular case the acquisition system is the same, the 
upper altitude PSF can be assumed to be carried out only 
by means of geometrical properties; other effects due to 
altitude increase (just like atmospheric effects, etc.) 
are not taken into account. In this hypothesis the time- 
scaling relation is used fRef. 3 pg. 32] i.e. if K is the 
altitude factor then 


(x) = hj^ (x/K) 

H, (f) = (Kf) 


( 2 ) 


If the two altitudes are double the other, the h^(x) can 

be represented as in Fig. 2.5b where = 1/2 f;}j; Fig. 

2.5c shows the convolution of this PSF‘ with ground object. 

Also the pixel dimensions are obtained by means of the 
scaling of the lower altitude pixel, i.e. now " 

2T, as it can be seen in Fig, 2.5a; the result of the sam 
pling is represented in Fig. 2.5e. It's obvious that now 
the samples' number is the half of the lower altitude one. 
The step corresponding to Fig. 2.5f and Fig. 2.5g is equi^ 
valent to that seen in Fig. 2.5f and 2.5g. 

It must be noted that in the frequency domain the sampling 
spacing is the same as for lower altitude (i.e. 1 /Tq) whi- 
le the samples number is reduced, as for space samples, to 
a half . 
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2,3. Determination o£ a scale -change function, 

A Fourier transform pair which permits one to carry 

out i 2 directly from ii must be defined. This pair must 
satisfy the relations 


.'32 

= i^ ® g 


C3a) 

C3b) 


Where with^lC^o) as denoted the DFT of ii(i 2 ) i-e* a com 
plex function, It immediately appears that such a determT 
nation can be easily made in frequency domain. In fact il 


Cu.v) = ^ i(u,v) 

exp [i 4^ 

Cu,v)j 

'J 2 (u,v) = I 2 (u,v) 

exp [ i 4> 

Cu,v)] 

(u,v) « G (u,v) 

exp 



Solving the equation C3a) gives ^(u,v) = G(u,v) exp ^i p(u,v^ 


where 


G Cu,v) = I, Cu,v)/Ij^Cu,v) 
p (u,v) = 4) (u,v) - tp (u,v) 


C5) 


► 


Liv and Gallangher [Ref. 4j give a detailed descriptio^n of 
the problems involved in minimizing the error between ^2 and 
'^l-h functions. The main problems concern ii and i 2 pha- 
ses ^choice; in the present form phases are, however, ne- 
glected so that the model is restricted to magnitude va- 
lues, i.e. 

G = G(u,v) = l 2 (u,v)/I^Cu,v) = i^amipi/j^dmp'L 


It is not so easy to solve equation (3b) if an operation 
inverse to the convolution were to be used. 

In order to realize a significative comparison, the same 
number of samples, for the two altitudes spectra, should 
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h?* used. Because has the primary component with 

N‘ samples, it's convenients to enlarge it up to N samples 
by adding zeros. Then the spectra to be compared are shown 
in Fig. 2.6. 



It must be noted that this enlargement has no consequen 
ces because the meaningiul part o£ the spectrum remain? 
unchanged; the spatial domain image is also enlarged but 
without any add of information. 

Referring to the primary component N samples, it's pos- 
sible to write 


jsampl 

2,enlarg 




jSampl 


true that 


OCf) * 

OCf) • Hj^(f) ' 
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jSampl 
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C8) 



where 
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HI ^2 Cf) - Hi ^2 

is the acquisition system continuous MTF sampled with 
1 /Tq step. 

The expression (8) states clearly that the G function 
can be carried out either directly from the data or from 
information on the acquisition system. 

Since spatial filters are passive devices, they don't am 
plify the signal. The G function, then must satisfy the“ 
constraint 


G C^iV) 1 

While this must be tested for different sensors, it is 
automatically true when H, is carried out by the scaling 
of H,*, this also pievents‘'any problem coming from Hi ze- 
ros . ^ 

The g function is obtained by means of the Inverse Discr£ 
te Fourier Transform CIDFT).'This g is constituted of N “ 
samples, but only few of them are non-zero, 

Before applying the proposed method for a generical scal^ 
ing-up, it is better to verify^it by performing^the fil~ 
tering in order to re-obtain 12 starting from ii, This 
operation can be performed either in the spatial or in 
the frequency domain. 

In the latter one, to be re-obtained, the multipH 

cation of G and must be realized. This gives, howe 

ver, if^eniarg necessary to reduce 

it by dropping (N-N')/2 zero value points at each edge. 

The filtering in the spatial domain is obtained convol 
ving 11 with the g samples. The output consists of N sam 

pies T spaced whose spectrum is just i^'^Q^iarg* ^ 

cimation of samples is made with a factor K equal to the 
scale factor (in the particular case of I 2.2 K = 2), 

jsarap just obtained. 

The decimation is in fact equivalent to the convolution 
in the frequency domain with a pulse train, spaced 1/KT. 


// 


This produces a loss of a part of the spectrum due to 
the aliasing; this part is however just the one containing 
zero values so the final result is similar to Fig, 2.5g. 

It must be osserved that the considered model is based 
on the assumption of the existence of a unique MTF for 
each altitude, i.e. in the hyphothesis of shift-invariant 
PSF. If this is not true, the model is no longer applica- 
ble in its present form. 
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3. SCALE EFFECT*. OPERATIVE APPROACH. 

The problems involved in the simulation of a higher alti^ 
tude image are: ” 

- determine the appropriate scale change filter 

- perform the filtering operation. 

3.1 Determination of the scale change filter. 

To solve the problem with a deterministic approach it is 
necessary no longer to consider data and to use only a 
priori informations. Referring to e.xpression (8) the G 
filter can be carried out by means of ratioing discrete 
MTFs. 

This obviously requires the exact knowledge of transfer 
function at each altitude; if these are not known, they 
must be constructed from the knowledge of imaging pro- 
cess, of all degradation sources involved and of their 
models . 

A quick survey of the principal degradations of an acquit 
sition system is reported in appendix II, but it appears 
evident that it's a very hard job to use this approach 
for carrying out the global transfer function )n the a^ 
sumption, however, that the lower altitude MTF is known 
and if the image seen at upper altitude bu the same acqui^ 
sition system is desired, it is possible to arrive easily 
at the G function. 

According, in fact, to (2), the G filter project procedure 
can be schemed as in Fig. 5.1. 


IS 
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It must be noted that the time -sea ling property can be 
indifferently applied either in space domain or in the 
frequency domain. 

The meaning of the operation G=H 2 /Hi is that explained 
in § 2.3. 

Neglecting the phase term, the function 


if ^2m,n ^ ^ 

''im.n 

CIO) 

if > 1 

^lm,n 


is computed. 

If the modulation transfer fmetion inf omiat ions are not avails 
ble, the scale change filter must be evaluated directly “• 
from the data. 

According to expression (6) the empirical approach is 
realised by ratioing the two altitudes images' spectra. 

In order for the use of the Fourier techniques to be mean 
ingful ,• ■ it is necessary to remove from the images all 
the shift-variant effects. Apart from the eventual SVPSF 
components of the acquisition system, there are some stan 
dard SV effects related to the acquisition itself. The 
actual procedure should consist in a deconvolution 
for these SVPSF, but empirical methods are usually used 
t'- correct data. 

In a first approximation two kinds of effects must be 
eliminated. The first is the so-called panoramic effect. 
This is due to the variation of the ground resolution 
when the view angle is changed (see Fig. 5.2). 


H. 


4ir. , n 


H 


lm,n 


m,n 
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I£ this geometrical effect is left uncorrected , it gi- 
ves incorrect spatial frequency values, especially for 
aircraft data where panoramic effect is more relevant . 
The correction can be done by using resampling techni- 
cpes which obviously yield an increase of samples number 
in the direction normal to flight line. 

Another effect depending on the position 

in the image plane is that due to the dif 
ferent atmospheric thickness for nadir and off-nadir ob 
servations. This yields a symmetrical lowering of the “ 
values at imngas edges respect to actual values. 

Instead of deconvolution approach, this effect can be 
statistically corrected by means of best-fit techniques. 
The assumption is made that for almost uniform regions 
in the image, the profile of the average values of the 
columns in this region can be fitted with a second order 
polynomial. 
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An example of this kind of correction is reported in 
Fig. 3.3 and in Fig. 3.4. The first shows a corrected 
and an uncorrected 3000m aircraft image, while Fig. 3.4 
shows the profiles of average values of the columns (a 
veraged on a number of lines) . Actually there are two 
competing atmospheric effects present in remote sensed 
data [Ref. sj . The first is the attenuation... which re- 
duces the magnitude of sensed radia-i^ion. The second is 
an active effect represented by scattering and/or emi^ 
sion. The balance between them depends on the direction 
of scan relative to the sun's position. If the scanner 
looks toward the sun, the dominant effect is the atte- 
nuating one. The opposite effect dominates if the scan 
ner looks away from the sun. it follows that the propo"- 
sed correction isn't completely exhaustive while a third- 
order curve should be necessary for a better correction; 
no attempt has however been made. 

Another effect is that coming from atmospheric thickness 
increase due to platform altitude increase. As it fol- 
lows from [Ref. 6j to simulate the increase of optical 
thickness, a linear transformation can be used which is 
independent from starting haze level, but dependent on 
added layer. It can then be considered shift-invariant 
and nothing has been done to take it into account. 

When all known shift-variant effects have been removed 
both from lower and upper altitude images, Fourier tech 
niques can be applied (in the assumption of SIPSF of the 
instrument) . 


In the following step Fourier transform of "corrected" 
images at different altitudes must be calculated. 

Given a discrete image in the form of the finite sequen 
ce {i(n^,n.,)} n^ = 0, . . . , N-l,n 2 ' = 0,...,N-1 its 

two-dimensional discrete Fourier transform is obtained 
[Ref. 7 ; Ref. 8] as the E-transform along the two circles 

•Wl -W2 

2 ^ = e^ , 2 , = e-' so that the transform is 


I = ^ |z^=exp[j C2mk^/N)j , S 2 =exp [j C27rk2/N)j 

N-1 .N-1 P 

= £ E i(n, ,nn)cxp I - j (2TTkTn,/N)i exp ; - j ( 2 TTk-n 2 /N): 

n,.=0 n.=0 ^ " L ^ J L J 
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with k^=0, . . . ,N-1, k 2 ® 0 ,...,N-l 

From an implementation point of view this two-dimensio 
nal transform is obtained by performing two one-dimen- 
sional transform as is shown in Fig. 3.5. 





Fig. 3.5 

It must be noted that the result has a phase which is 
such that the meaningful parts of the transform are at 
the corners of the display. To help visual examination 
of the transform, it's convenient to shift the origin in 
the center of the N x N points frequency domain. This 
is obtained either by re-ordering transform coefficients 
or by making use of the translation property of Fourier 
transform pair 


f(x,y) exp [ j27rCuQX+VQy)/Nj<J=(^ FCu-u^,v-v^) 


i.e . 


( 12 ) 


i(nj^,n 2 ) exp [ j (2Trk^n^/N)] 

If the shift is k^=k2=N/2 
exp f j2ir (k^n^+k 2 n 2 )/Nj = 


exp f j 2 iTk 2 n 2 /N)j<j=(>^(mj^-k 2 ^,m 2 -k 2 ) 
then 

gjTr(nx+n2) = (-.iynx+n2 

(15) 


// 


fig. 3.3a 



uncorrected 


fig. 3.3b 
corrected 




fig. 3.4 
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so that multiplying pixel of cooordinates (n^,n 2 ) by 
(_ 1 )‘^ 1'*’*^2 j the requested transform is obtained. 

The software in its present form (implemented on IBM 
370/1S8) produces only the spectrum of the image, i.e. 
being generally the Fourier transform complex 


F(u,v) » R(u,v) + j I(u,v3 (14) 


the available data are 


F(u,v) 


[r2Cu,v) . I^(u,v)] 


CIS) 


After the two altitudes images' spectra are obtained, in 
accordance to (6) the scale change transfer function is 
carried our by ratioing these spectra. This ratio can be 
performed either on a point-to-point basis or in an in 
tegrated manner. 

The first approach is certainly more general, but it is 
also more subjected to local fluctuations. The integra- 
ted approach consists in subdividing the original spec- 
tra in a lot of little areas as it is shown in Fig. 3.6. 
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The ratio of partial energies of corresponding regions 
gives the shape of the G function. On the assumption of 
circular symmetry of the scale change filter, the regions 
to be ratioed can simply be the rings, 

With respect to the proposal to obtain the scale change 
function for whatever altitude increase, two cases can 
occur: 

i) If only two altitudes data are available, it is pos- 
sible to determine only one G function. Supposing the 
acquisition system used is the same, in order to si- 
mulate other altitudes data, the system MTF must be 
known. It's then necessary to proceed by attempts with 
the following steps; 

- the hypothesis is made that the system MTF belongs 
to a certain class of filters; 

- by means of time-scaling and of the model of § 2 

the corresponding scale change filter is carried out; 

- a comparison between supposed and true scale change 
function is realised which gives an answer to how 
good the choice is. 

When the system MTF is carried out, other altitude in 
creases can be simulated with the method shown in FigT 
• 5 . 1 * 

Just as an example two filter classes (Lanccos-Cap- 
pellini and Butterworth) are presented in Appendix III. 

ii) If at least three altitudes data are available, it is 
possible to determine at least two G functions. Writing 
these function using the altitude factor as parameter, 
it is possible to obtain directly the G filter for wha 
tever altitude increment. 

3.2. Filtering implementation. 

As it is well known, the filtering operation can be per- 
formed either in the spatial or in the frequency domain 
the choice depending on image site and impulse response 
site. If non- zero impulse response array is large, Fourier 
domain is more computationally efficient while direct pro 
cessing is advantageous for large site image. For very lar 
ge images a block mode filtering technique can be used 
which processes separately adjacent overlapped blocks 
fRef. 10 chapt. 11 j . 


Si! 


- 19 - 


i) Spatial domain 

The two steps of direct convolution and points' deci- 
mation can be grouped into a single operation of "con- 
volution by step" on the basis of reduction factor. This 
permits only integer factor reductions so if it is not 
case, it is necessary to apply resampling techniques; 
this is however equivalent to use, as scale change fun- 
ction, a filter ^different from the projected one. 

Attention must be used in handling edges data in order 
to have a meaningful convolution also along the borders, • 

ii) Frequency domain 

In this case the filtering operation consists in multi 
plying the image transform by the corresponding samples 
of scale change filter Cmultiplication of modules; dif- 
ference of phases) . Now there is no more the restriction 
to integer reduction factor. While, in fact, in the spa 
tial domain decimation must be used, in frequency domain 
the same result is obtained cutting away an edge of spec 
trum data so that, depending on the size of this edge, an 
arbitrarly reduction can be obtained. 

4. COMPARISON TECHNIQUES BETWEEN TRUE AND SIMULATED ‘ 
IMAGES . 

After the filtering has been performed some methods must 
be established to carry out the comparison between simu 
lation output and true image. A pixel- to-pixel comparison 
must be disregarded because such an approach presupposes 
a very accurate registration between images. Also taking 
into account the phase factor and on the hyphotesis that 
images don't have any geometrical diitorsions, the re- 
gistration averaged error for N reduction factor is 1/2N 
of the low resolution pixel diinensi’cn. It follows that for small 
reductions (just like aircraft data simulations) this er 
ror would be rather high, 

Instead of local techniques it is then convenient to use 
the image global informations. These techniques can be 
summarised as follows: 

i) image analysis in the space domain 

- histograms comparison 

- entropy computation for the two images. Referring 
to rRef. lit it is possible to define the entropy 
of an image as 





( 16 ) 


N 

S - - ^2^ p.log2Pi 


Where with p^ has been denoted the level probabili- 
ty (it is interesting to note the introduction, in the 
same paper, of a specific entropy, i.e. a value indipen- 
dent of 'image area) . _ 

This method permits to realize a comparison of two images 
simply by the comparison of two numbers. 

ii) Image analysis in the frequency domain. 

- computation of spectra correlation. 

In the transform domain the problems arising from 
a point-to-point comparison are smaller than in the 
space domain. Here, in fact, the registration of two 
sijectra gives no problem unless the original lines 
of flight were not perfectly parallel; in such a ca- 
se the true and simulated image spectra are slightly 
rotated one in respect to the other. If this is not 
the case, it is then possible to estimate the simi- 
larity between Images by means of a measure of cor- 
relation of spectra. Anyone of known fidelity fun- 
ction [Ref. 12, vol. 3.1, 1.4.2.3.J can be used; the 
sum of absolute values of differences will be used 
in the following 

- Spectra energy computation in corresponding regions. 

If an accurate registration between spectra is not 
possible, an averaging technique must be used. Remem 
bering that the energy of a signal is the square of 

its spectrum, the total energy of the Fourier transform 
of each image is calculated adding the energy at each 
point (u,v) 


N-1 N-1 

E™ = S S ECu,v) (17) 

^ u»0 v-0 


with u,v=0, . . . ,N-1 . 

Calculating then the energy enclosed in circles or 
circular rings with origin at the center of the fre- 
quency domain, it is possible to compare the corre- 
sponding energy percents S of the two spectra . 

This is obtained by calculating 


e - 100 


r 


1 


(18) 




Z Z E(u,v)/E,r.! 
u V 


where the summation indices range over the samples 
lying inside the circles or the rings. 

Another averaging technique makes use of circular^ 
sectors. If, however, a rotation between spectra is 
present, the circles/rings method must be preferred. 

5. APPLICATION EXAMPLE 

Just as an example of the use of techniques seen in pre>^ 
vious paragraphs, a simulation of an aircraft day image in 
the infrared channel is now examined. The data collection 
was made during the 1973 joint flight experiment on test 
site Sibari by means of a Bendix MSS. During this compaign 
data were detected at different aircraft altitudes (500, 
1000, 1800, 3000m) but as a first attempt only 500 1000m 

simulation had been put into practice. 

The region considered (seen at 500m) is reported in Fig. 

5.1. 

According to the procedure described in § 3.1, the following 
steps have been implemented: 

- correction of 500m and 1000m original images. 

Since the Bendix PSF isn't known, the only corrections 
have been made concern the optical depth increase for 
scan angle and the panoramic effect. About the former, 
the correction has been at first calculated referring 
to an almost uniform region and has then been applied 
to the whole 500m image, Fig. 5.2 shows the profiles 
(averaged on a number of lines) of the uncorrected and 
Corrected image. The same procedure has been repeated 
for 1000m image. The panoramic effect has been corrected 
resampling data (with a cubic polynomial) in the scann- 
ing direction. The algorithm has produced, using 803 
input samples, 1256 output pixels. A visual effect of 
this correction is shown in Fig. 5.3 where uncorrected 
image (Fig. 5.3a) is compared with the corresponding pa 
noramic effect corrected image (Fig. 5.3b). The two ima 
ges are left-aligned (corresponding to flight line) so 
that pixels' number increase becomes evident (a symme- 
tric effect on the other side) . 




IhO 00 .20 CC 








- 22 - 


- Choice o£ test area. 

The simulation procedure is performed on a segment 
512x512 pixels dimensioned. 

In the whole region of Fig. 5.1, the area shown in 
Fig. 5.4 has been chosen because of the presence of 
spatial structures. For homogeneous area, in fact, the 
scaling-up output is independent from used PSF so this 
case isn't meaningful for test. 

The scale change procedure should give an output simi- 
lar as much as ia possible to the true 1000m image 
ported in Fig. ^.5). 

- Fourier spectra computation. 

The developed software in its present form, neglects 
the phase information; besides since the dynamic range 
of transform is very large and in order to display it 
in image form on 1 2s imaging system, it's necessary to 
Compress the coefficient values in byte format. The 
whole spectrum range is then linearly subdivided into 
256 levels. 

Images of Fig. 5.4 and 5.5 are used as input to Fourier 
spectrum computation program. The resulting outputs are 
respectively reported in Fig. 5.6 and 5.7. For visual 
reasons they report the quantity 1+logio |Hl* 

As it has been stated in § 2.2, the spacing between fr^ 
quency samples is the same, i.e., Au-1/Tq where T0=NAx7 
In this case for 500m image spectrum ax- 1.2Sm, N=512 
while for 1000m image spectrum, ax- 2.5m, N=256. It 
follows that the frequency domain dimensions are (see 
also [ Ref. 8 pg. 7 8] } 

512 du 

256 du ■ 


The maximum revealable frequencies are 

1 

0.5 0-*^ cycles/m 

0-5 — .= 0.2 eye 1 e s / m 
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Ratio of spectra. 

As it has been stated in § 4, the ratio of spectra can 
be realized either in an integral mode or in a point- 
to-point mode. 

An example of the former approach is reported in Fig. 

5.8 which shows one of the rings the Fourier plane can 
be subdivided in. The technique used is, however, the 
local ratio one; the (normalized) output is in Fig. 5.9. 

The presence of structures in the spectrum can be ascri. 
bed to a little non-compensated rotation between the 
two original images spectra. 

Determination of an analytical scaling-up filter. 

As it has been stated in §3.1, it is convenient to de- 
termine which class of analytical filters fits the true 
scaling-up filter. 

In this case it has been seen that the ratio of a But- 
terworth filter of order 3 and cut-off frequency 120 
with its scaled version produces a G' function which r^ 
sembles the shape of the G filter (Fig. 5.10). Accord- 
ing to the model, the acquisition chain object-data at 
500m altitude can be then represented by a MTF shaped 
just like a Butterworth filter order 5) so the 

corresponding G filter for whichever other altitude in 
crease can be easily carried out. 

Scaling-up simulation. 

The image filtering has been performed in the space do 
main. Frequency domain filtering requests, in fact, the 
multiplication of the spectrum by the filter. Since both 
have a quantization error and because the result must 
be re-quantized and antitransformed, referring to the pre 
sent software, this approach would be meaningless (for 
a good treatment of quantization techniques see [Ref. 9 
chapt. 7J ) .The two dimensional spatial array obtained 
from G' antitrasformation has been convolved with the 
image of Fig. 5.4 and the result has been compared with 
Fig. 5. 5. 

True and simulated images comparison. 

Referring to the techniques exposed in § 4 , the histograms 
of the whole images have been obtained (Fig. 5.11). As it 
is evident, there is no similarity between them. This e_f 
feet becomes, however less evident if the histograms re 
garding a region symmetrical respect to the flight lines 
are considered (Fig. 5.12). The difference still remain- 
ing between them can be further reduced by means of a 
linear transformation. Calling, in fact, xi,Mi,S]^(X2,M2>^2^ 
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respectively the gray level, the average value and the 
variance of the simulated (original) image, a new ima- 
ge has been obtained with gray values carried out 
by means of the relation 


where 


X 


aXi + b 


a ® 


S, 

b = M, - aM, 


Fig. 5.13 shows the comparison between the original and 
the new simulated image histograms. 

The hypothesis can be made that tiiis linear transforma- 
tion accounts for atmosphere increase with altitude, 

When the transformation has been applied to the whole 
image no improvement has been obtained. 

Two possible causes are suggested to justify this re- 
sult: 

- along-flight direction pixels correlation (the cor- 
rection of panoramic effect has been realised, in 
fact only in the scanning direction without conside- 
ring the enlargement of the pixel also in the flight 
direction) 

- shift-variant components present in the Bendix MSS 
PSF. 

As a further comparison method, the simulated image 
(complete) Fourier snertrum has been compared with the 
original one. From a visual inspection the 

filtering output is too lacking in high frequencies in 
respect to the true 1000m image, i.e. the filter used 
has* been revealed as the incorrect one. 

As an experiment all the operations have been repeated 
with the Lancsos-Cappellini filter with parameters 
{N=16; a-0.25; m=15} . Histograms and spectra reveal 
however, no meaningful improvement respect to the first 
attempt. 


3 / 


§‘ I » > 
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6. CONCLUSIONS 

The most important results arising from this work can 
be summarized as follows: 

- singling out of main problems which must be faced in 
scaling-up simulations 

- determination of a simulation method 

- software implemen-tation of the method. 

The lack of informations on used acquisition system 
hasn't permitted the realization of an exhaustive test 
of the proposed scaling-up methodology. 

Even more so HCMM images haven't been not simulated also 
because of the difficulty of registering the few pixels 
coming from test-site data with a whole HCMM scene. 

Further developments can be however proposed to improve 
the method. They can be so summarized as follows: 

- utilization of atmospheric models to take into account 
the increase of the atmospheric layer 

- correction of panoramic effect also in the flight di- 
rection 

- better registration between the two altitudes starting 
images (if a rotation must be performed, it is convenient 
to subdivide it into two opposite half-rotations; this 
avoids the resampling affects in different manners the 
two images) 

- introduction of the phase term in the Fourier transform 
computation and usage of non-integer values for the 
transform itself . 

Beside these, the module could be developed to carry out 
the scale change filter in a parametric form (see i 3.1) 
using three altitudes data; this avoids, in fact, the step 
of the determination of instrument MTF. 


Many aspects of the proposed empiric approach make useful 
the* approach itself only for the simulation of different 
altitudes data detected by the same acquisition system. If, 
on the contrary, images detected by different systems must 
be simulated, only the deterministic approach can be used. 
The problem must then to be faced to know with a sufficient 
approximation degree the MTF of the acquisition system and 
after to apply the method of pg.l2 
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APPENDIX I 


SURVEY ON POINT SPREAD FUNCTION (PSF) AND MODULATION TRAN 
SFER FUNCTION (MTF) . 


The image formation process by means of an imaging system 
(see Fig. A.I.l). 




Fig. A.I.l 


involves a redistribution of radiance coming from the object. 
Generally this redistribution must take into account the ra 
diance in the whole image plane, but if the process is linear, 
then the image of an extended source is obtainable by super 
imposing the images of single object points. The degradation 
introduced by imaging system is then represented by the spa- 
tial distribution on image plane of the single source point 
and it is mathematically expressed by the point spread fun- 
ction (PSF) . Denoting with o(x,y) and i(x,y) respectively 
the object and image functions and with h(x,y) the PSF, the 
linear imaging process is then represented by 


f -f* co ^4*00 


i(x,y) = 


— CJO/ 


- 03 


o(C,n) h(x,y;5,n) d^dn 


(A.I.l) 


The four indices for impulse response denote that its 
shape can change depending on input and output coordi- 
nates. If, on the contrary, each point in the object pla 
ne contributes to the image point with a "strenght" de^ 
pending only on its relative position, the system is said 
to have a shift-invariant point spread function (SIPSF) . 
In this case the (A.I.l) becomes the superposition inte- 
gral 


i(x,y) 


^ 4*00 ^ 4*00 

J -ooJ -00 


o(C,n3 


hCx-?,y-n) 


d§dn 


(A. 1.2) 


This expression is nothing but the convolution of o(x,y) 
with h(x,y), denoted also with 


i(.x,y) « o(x,y) ®e h(x,y) (A. 1. 3) 


If the Fourier domain is now considered (obviously with spa 
tial frequencies), it is possible to associate each spatial 
function with its Fourier transform. These new functions are 
complex ones and so each of them can be written in the form 

R(u,v) » P(u,v) + jQ(u,v) “ M(u,v) (A.I.l) 


Starting from the Fourier convolution relationship (i.e. the 
corrispondence between convolution and multiplication in spa 
tial and frequency domains as it is shown in Fig. A. 1.2) ~ 
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Fig. A. I. 2 

the expression (A. I. 3) becomes 

I(u,v3 = 0(u,v) H(u,v3 (A.I.S) 


Considering the PSF transform, it is usually known as Optical 
Transfer Function (OTF) . In particular its amplitude respon- 
se is Called Modulation Transfer Function (MTF) while the pha 
se response is the Phase Transfer Function (PTF) . 

It must be noted that in the practice phase variations are 
often disregarded and so the imaging system degradations 
are represented either by PSF in the spatial domain or by 
MTF (referred to in the text by means of capital letters) 
in frequency domain. This pair of functions is usually cal 
led Fourier pair and it is sufficient to identify the imag- 
ing system. 




A useful characteristic of MIT is its "multiplicative 
property", i.e. if an imaging system is set up with more 
cascaded devices, the global MTF is carried out by the 
product os single components’ MTFs. 

It must be underlined that the definition of an univocal 
MTF is possible only for SIPSF. 

If on the contrary shift variant PSF (SVPSF) are involved, 
Fourier techniques are no longer applicable and the use 
of very complex mathematical techniques is required [Ref. 
25 pg. 24, pg. 52 ] . 
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APPENDIX II 

SURVEY OF PRINCIPAL DEGRADATIONS OF AN ACQUISITION 
SYSTEM. 


Referring to a typical remote sensor configuration, 
it can be represented as in fig. A.II.l 



If its performance is to be evaluated, two approaches 
are possible: 

- the direct method, i.e. the Modulation Transfer Func 
tion of sensor system is determined by a variety of 
laboratory measuring techniques [ref. 13 pg. e-Z97j 

- the indirect method i.e. a model of the system is 
constructed Por this approach the mathematical model 
of each component of sensor system must be known and 
then the global MTF is carried out by making use of 
multiplicative property for MTFs of cascading system 
elements [tef. 13, pg. e-29l] 










For the indirect method to give meaningful results, 
it is necessary to know all the degradation sources 
of the system and their mathematical formulation. 

Examining at first the sensor system, its single com 
ponents can be considered: 

- scanning system 

Quite apart from the scanning mechanism used (rotat- 
ing or oscillating mirror) there are often variations 
in the angular rate of the mirror which distort the 
image changing its frequency content. 

It must also be remembered that the ground is scanned 
an_ angular displacement to either side of the nadir. 
This gives a change of effective pixel's dimensions* 

A good account of this problem can be found in ref. 14 
pg. 33J . It can be summarized here saying that in 

considering only frequencies in the scan direction 
(fy=0), the corresponding transfar function is 


TapCf^-) * sine (D^f^ sec-^ ) sine (Dyfy tan^seci?) 
where D = W H/F 

A , / A , / 

W" = detector dimensions 

A , / 

H = platform altitude 

F = focal lenght 

t 9' = off-nadir angle 

In a general manner, called ~ ^ ^y sec^i^J 

a generical spatial frequency measured in cycles/radian , 
it is possible to write 




where 6 are the angular dimensions of IFOV 
- optical system 

Since the optical system has finite dimensions, the 


(il Bril')*# M 
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diffraction phenomenon is unavoidable. If the use 
conditions., allow the assumption . that the system 
is aberration free and if the aperture is circular 
shaped, the PSF.has the usual Airy diffraction pat- 
tern with the central disk (i.e. limited by the first 
obscure ring) containing the 84‘i of global radiance. 

If the resolution properties are expressed as. trans- 
fer function then [ref. 14} 


T 


optics 




^ arcos a - a (1-a 0 c a < 1 

0 a 5 1 


where a = f„ X/D 

D = lens diameter 
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It must however be remembered that the usual scanner 
configuration has the primary mirror partially obscured 
by the secondary mirror and this gives an energy redi- 
stribution in the diffraction figure. The theory of 
this can be found in {’ref. IS; ref. 16, pg. 336j ; 

to have, however, an idea of the behaviour, if p 
and HP are respectively the external and internal 
radii of the aperture, fig, k.ll.l shows the dependence 
of MTF from n parameter. From this it is evident j that 
the introduction of a central obscuration reduces the 
responce at lower frequencies while it improves at 
higher frequencies. 





fig. A. II. 2 

Up until now aberrations have been neglected; this is 
however not quite correct because in the components 
of photoelectronic imaging system they are the major 
degradation sources. 

Referring to aberration as a phase variation, it is 
possible to demonstrate that it produces at each fr£ 
quency a lowering of transfer function, in a general 
manner, in the phase variations expressions polynomials 
appear whose degree characterizes the aberration 
itself. 

The linear terms are usually known as distortions. 
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Their treatment . is made rather complex by simulta 
neous dependence on the cube o£ the field angle, “ 


The second degree term represents a focusing error . A 
good account of this degrading source can be found in 
[ref. 17] ; the relative transfer function for circul 
ar aperture has the expression [ref. 14] ” 


T 


op 



tt/z 

’ Li sin U + cos U - 1 j.,. 

— y— - •d<|> 

■'6 (C cos )‘^ 


where 

U = C j^(l - sin^ ij; 
C = ir (D/F)^ a e/X 


a cos 'll 


cos i|; 


a » X/d 
D = lens diameter 


F = focal lenght 

e = displacement of the detector from the focal 
plane 
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Its dependence on the e parameter is shown in fig* A^II.3* 


errors 
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fig. A, I I. 3 


For values e > 2 X , the function has the sine trend 
The cubic term aberration is the so-called "coma". An 
example is reported in fig. A. II. 4 




\r:Ol 





Fig. A. 1 1. 4 





It must be remembered that the related MTF is no 
more rotationally invariant as is evident from fig. 
A. II. 5 
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fig. A. II. 5 


which represents the MTF for a IX coma; the equal 
modulation lines are 0.05 spaced. For an interesting 
account of these two degradations see [vef. 18J . 

The fourth degree term is the spheric aberration. 

For their intrinsic nature, the aberrations have 
space variant PSFs; it follows that for a complete 
OTF representation, it would be necessary to give 
separate plots depending on view angle, focal posi- 
tion, wavelength and so on. 

It is interesting to note that by using polar coordi- 
nates and with a suitable change of variable in r, 
it is possible to reduce some aberrations ii.e. coma) 
to linear spatially- invariant degradations fref . 13, 
Pg- 18] . 

- Field-defining system 

Under the same use conditions the aperture shape in- 
fluences the system response as it appears evident 
from fig. A. II. 6. - 
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fig. A.n.6 

- Detecting syiitem and electronics 

It must be remembered that electronic system converts 
a spatial function into a temporal one. Referring only 
to the case the relation is Tref. 14j 

f = f H 'S^ sec“i9' 
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it follows that detector and electronic characteristics, 
expressible in terms of temporal frequency transfer 
functions, influence also spatial frequency MTF. 

Since the detecting system can be based upon different 
physical processes (photo-multiplyers , photodiodes, 
photoconductors, etc.) also MTF varies depending on 
the model; an example of MTF derivation is given in 
[ref. 19 pg. S3] . 

- Atmosphere 

In respect to the image which would be detected by the 
sensor in absence of any medium, the presence of the 
atmosphere produces a blur The effects of this degra- 
dation source are of two kinds: 

- a radiometric effect coming from scattering, absorp- 
tion and emission of radiation; 

- an effective geometrical blur, representable in the 
form af a transfer function, caused by atmospheric 
turbulence . 

Examining the latter degradation, it is clear that 
there is a significant difference for the case of 
short and long time exposures For the long one, it 
has been demonstrated that the blurring filter may be 
approximated by the function H(u,v) = exp f -G (u2+v2) 5/6J 
where C is a constant depending upon' the atmosphere 
turbolence . 

For convenience, the term 5/6 is often replaced by 
unity so a gaussian curve can be used. A detailed 
account can be found in the original paper by Hufnagel 
and Stanley 20] . 

As it appears also from [ref. 14] , in the complete 
MTF expression there is also an empirical form inte- 
gral whose numerical value depends from platform 
altitude . 

Among the other degradation sources not still mentioned 
there are, for example, vibrations and, generally, the 
motion of the platform. This topic isn't considered 
here buffer, a good account of transfer functions coming 
from image motion see [ref-. 14 pg. 37; ref, 16, pg. 436; 
ref. 21; ref. 22, pg. 101] , 
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APPENDIX III 


LOWPASS PILTBRS; SOME EXAMPLES 

As ic is known, the lowpass filters are characterized 
by attenuating more the high frequency content of an 
image than the low frequency one* Eor the usage these 
filters are provided for in this work, only filters 
belonging to the FIR (finite impulse response) class 
are considered. This means that the impulse response 
is of finite extent with zero values outside some fi- 
nite limits. 

Furthermore, considering the system is memoryless, the 
filter output samples must be determined only by input 
samples i.e. nonrecursive filters are involved. The 
filters result so to be zero-phase-shift-filters, 
i.e, they don't alter the phase of the transform. 

The unit sample response for many types of two-dimen- 
sional filters is of infinite extent; to satisfy 
then the FIR request, a truncation is needed. Using 
as truncation function t!‘ ' rectangular one, this re- 
sults in the introduction in the frequency domain of 
a percentage of ripple and overshoot. It can be seen 
that using window functions having in the non-zero 
region values different from unity, it is possible 
not only to reduce the error, but also to suitably 
shape the MTF. A lot of window function can be carried 
out to proiect filters as it comes out from fref. 23 
vol. 2 

Just as an example tlie Lanczos-Cappellini window is 
more closely examined [ref, ?, pg. 6§J . 

The starting point is the ideal 2-D circularly symme- 
tric low-pass filter with frequency response 
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where wq « In fj^.) i.s the cut frequency shown in 
fig. A, II 1 . 1 . 
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fig. A.III. 1 

The corresponding discrete impulse response samples 
are given by 

J, a »4?+i4) 

h(n,,nO = a — - — =- (.A. III. 2) 

X i. r-:^ 7 

/nj + n; 


where Jj^: first order Bessel function 

a = normalised frequency fj./fg with fs sampling 
frequency. 


This function is clearly not limited and with a notice 
able ripple. The window proposed by Lancsos in the one 
dimensional case and extended to two-dimensional one 
by Cappellini is 



where m is a real positive number and 


n 


1 ’ 


n. 


N-1 



The final impulse response is then carried out by 






multiplying expressions (A. III. 2) and (A. III. 3). 

It is interesting to note that changing the set [N; 
a ; mj it is possible to act on MTF characteristic. 

An example is given in fig, A. I I I. 2 and fig. A. I I I. 3 
where the MTF dependence from respectively, mi and N 
is shown. 

The window method is usually implemented in the space 
domain; it presents however some difficulties as, 
for example, in determining the impulse transform or 
to locate the band edges. 


An e.xample of filter carried out in frequency domain 
is the Butterworth lowpass filter (BLPF) . 

The BLPF of order n and with cut off frequency f 
is defined as 


HCf) * i (A. III. 4) 

1. [f/f J 2n 


where in the two-dimensional meaning f ' ‘ 

\ y 

Here the cut-off frequency has no longer the meaning 
of fig. A.III.l, but it becomes the frequency for 
which H(f) value corresponds to a certain fraction 
of its maximum (usually 0.5 or 1 // 2 ). 

One of the main advantages of this filter respect to 
the former is the possibility of controlling the degr 
of attenuation in the transition band. 

To carry out the discrete version of this filter, 
the frec}uency sampling method could be used [ref. 25 
§5.2.2.] or another of the techniques described in 
ref. 24"” §2-5] ; to obtain, however, a filter resembl- 
ing the BLPF in a manner as simple as possible, the 
expression (A. II I. 5) can be used; 
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fig . A . Ill . 5 . 
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V 


H(n2,n2) 



(A. III. 5) 


where with Af frequency sampling inter- 

val . 


The cut-off frequency makes the H value one half of 
its maximum value. 

Fig. A. III. 4 shows the trend of this class of filters 
varying the cut-off frequency while fig. A. III. 5 shows 
the dependence on filter order. 

Further examples of filter can be found in fref. 8 
§ 4 . 5 . 2 ] . 
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